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BESSEL FUNCTION - Jn(x) - SUBROUTINE (FLOATING POINTS) 

1. Puipos e 

To provide a subroutine to compute for a given argument x, all 
significant orders of Bessel functions of the first kind Jn(x). 

2. Restrictions 

The argument x must be in the range 0<rx <207«0631» 

3. Method 

3,1 From the value of the argument x, a starting order is determined 
from the relationship: 

1 
ng = 7 + (65.^lBx + 1.1765x2) (truncated to a multiple of h) 

The value ng, which is given in Table I as a function of x, 
is the total number of orders generated for argument x, and may be defined 
as the smallest multiple of h that satisfies the inequality: 

J,^(x) < 1 X lO'-'--'- for all n ^. ng 

Arbitrary values of 2 and 0, respectively, are assigned to 
JngCx) and J^ +i(x). The subroutine then computes all orders down to 
and including Jo(x) successively from the recurrence relation: 



Jn-l(x) = (2ll) Jn(x) - Jn+i(x) 
\ X I 



The resulting table of ng unnormalized Bessel functions, from 
Jq(x) to Jns-l(x), inclusive, is then corrected by dividing each entry 
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by a normalization constant N determined from the unnormalized values in 

the following way: 

N = Jq(x) + 2 ^ J2n(x) 
i=l 

with the summation terminating with the term Jn^Cx). Actually, in the 
subroutine, the summation is accumulated as the unnormalized Bessel func- 
tions are generated, 

3.2 Accuracy 

The accuracy of the Jn(x) functions obtained with this subroutine 
depends upon the portion of the particular Bessel function in question, 

3.2.1 ■ Monotonic Region (x<n) 

Since in this region the J^{x) functions contain no zeros 
(except at the origin), the concept of relative error is to be preferred. 
The relative error in this region is smallest at the order nearest x, and 
gradually rises smoothly, for successively higher orders, until an order 
about $ or 6 below the highest computed order (ng-1) is reached. From 
this point to the highest computed order, the relative error rapidly 
increases and is largest for the highest order, 

lypical errors for this region are as follows: 

RELATIVE ERROR FOR VARIOUS ORDERS (x<n) 
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-"■Last place agreement within two units. 
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Lack of adequate reference tables prevents an extension of the 

foregoing table to higher values of x, 

3.2.2 Oscillating region (x>n) 

•In this region, the Jn(x) functions are periodic with distances 

between zeros of 3-9 units in x depending upon the order. Thus absolute 

error appears to be a more useful measure in this region. 

The absolute error of the J^Cx) functions as computed by this 

subroutine is highest in the vicinity of the zeros, and lowest in the 

vicinity of the maximums and minimums. An analysis of the errors in the 

vicinities of the zeros reveals that the average error obtained may be 

expressed 

laverage absolute error | -^ (3 + 0,21x)lO"' in the vicinity 

of the zeros 

The errors obtained in the vicinity of zeros have considerable 
scatter, and maximum errors may be twice the mean error approximated above. 
In addition, when the derivatives of the function in the vicinity of the 
zero is positive, the calculated values given by this subroutine are nearly 
always greater than the true values and vice versa. Although the above 
expression for the zero-vicinity error was derived principally from 
Jo(x) and JpnCx) data, it appears to be applicable for all orders given 
by this subroutine. 

The average error in the vicinity of the maxima and minima of 
the J (x) functions appears to be more or less independent of both x and 
the order n. Typical absolute errors have been observed to be about 
3 X 10" with maximum errors less than 6 x 10"-'"'-. 
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h» Usage 

i;.l The argument must be a floating point number in the A and R 
registers • 

U.2 ' Calling sequence 

tn. ) FCA argument 

TRA Lq ■ 
o< + 1) PZE NMAXg 
PZE KREQ_ 



C!i + 2) Error return 
Normal return 
ii,3 Explanation of symbols 

U.3.1 MIAXq is an octal interger (at binary point l8) that serves 
to specify the storage allocation of the J (x) table. MAXq normally has 
the value l^OOg, which provides for storing the maximum of 256-j_q orders 
(from Jq to ^2^5^ inclusive) designed for this subroutine, Jq is always 
stored at location Lq + 0300q, J-j at location L^ + 0302q,.,.,J2^^ at 
location Lq + 1276qo 

If a user knows that his arguments will never exceed 100, for 
example, then reference to Table I shows that l^O^^o orders will be gen- 
erated (from J to Ji39-,Q, inclusive), or in octal 2li|n orders will be 
generated. Thus by setting MAXs '■= 2liig, twice this number (2 locations 
are required per order) or lj.308 locations are reserved for the J^ table. 
Two objectives are therefore accomplished: 

a. Memory locations starting with location Lq + 07308 are 
unused and are thus available for other purposes. 
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b. Arguments in excess of 101. [;9ii^ (from Table I will 
cause an error return, 

.!l,3.2 This subroutine normally computes all significant orders 
for a £;iven argument x, with the total number of orders computed given 
in Table I as a function of x. For this normal type of operation, NREQg, 
an octal integer (at binary point 38) should be set to zero. 

But if, for example, a type of problem arises where say the 
first 20-,Q orders are always required, regardless of the value of x, then 
NREQ may be set to 002l4g. The subroutine calculates the normal starting 
order, and then raises it, if necessary, in jumps of h to provide a 
starting order of from h-l higher orders than that specified by NREQ. 
This assures that all orders up to and including that specified by NREQ 
will be obtained to the full accuracy of the subroutine, 

Iul|#l Error return is implemented when the argument is zero, 
negative, or exceeds the limits imposed by NMAX (see example in l;,3.1.b). 
The error return is made with the argument x in the A and R registers, 

k>h,2 Normal return is made with the address of the highest order 
computed in the accumulator at binary point 38. The highest order itself 
(ng-l) is displayed as a decimal integer, 

h»S The subroutine is relocatable. Relocation subroutine RUG 1075 
is included on the tape, 

ii.6 The subroutine destroys the L and V loops, 

ii,7 Options available are discussed in i;,3.1 and ii,3.2, 
^« Codin g Information 

5.1 The subroutine occupies three tracks or 192]^o words. The J^ 
table always starts with J at location Lq + 03008 and extends up to a 
raaxiraum value of J25^ at Lq + 12765, 
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5*2 Unused locations. A total of 13 locations are unused as 

follows: 

L + 0C03 
° 0006 
0021 
0022 
0101 
0111 
0121 

0127 
OlUl 

0151 

0215 
0223 
0277 

5.3 Locations used for temporary storage total 10 as follows: 

Lo + 0010 
0011 
0013 
0023 
0075 
0076 
0077 
0150 
0251 
0252 

5.i4 Constants (floating point) 

1.1765 L + 0136 

2 ° 0216 

65.518 0226 

7 0270 

Constants (fixed point) 

10 at 39 L + 0212 
021ii 

^^^ Execution time. The execution time for the subroutine depends 

upon the number of orders generated and is approximately as follows: 

t = ng/lO where ng is the number of 
orders generated 
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6, Checkout 

To check out this subroutine, numerous comparisons were made 
with the classic tables of Bessel functions. The accuracy figures given 
in 3,2 stemjned from an analysis of these comparisons. Reference tables 
include: 

a. Meissel's l8-decimal table of J^{x) at unit interval in 

X up to X = 2ii as reprinted in Gray, Mathews, and I^cRobert, 
"A Treatise on Bessel Functions," second edition, London, 
Macmillan, 1922, 

b. Harvard Computation Laboratory. Tables of the Bessel 
functions of the first kind, 12 vol. 19ii7-19^1. Gives, 
for X less than 100, l8 decimal values for n«0ton=3, 
inclusive, and 10 decimal values for all higher orders. 
Interval in x ranges from 0.001 to 0.1. 
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TABL3 I 



0. 

0.3790 

1 .2100 

2.H699 

H.IO80 

6.0695 

8.3018 

10.7580 

13.3980 

16.1889 

•19.1037 
22.1207 

'25.2226 
28.3952 
31.6273 

•34.9097 
38.23^^8 

k] .^^^66 

i|8.i+io6 

51.8551 

55.3205 

58.8042 

62.3042 

65.8186 

69.3458 



8 
12 
16 
20 
24 
28 

32 
36 
40 
44 
48 
•52 

56 
60 
64 
68 
72 
76 
80 
84 
88 
92 
96 
100 
104 



69.3458 
72.8845 
76.4336 



79.9920 

83.5587 
87.1332 
90.7145 
94.3022 

97.8957 
101.4945 

105.0982 

108.7063 

112.3186 

115.9348 

119.5545 

123.177^ 
126.8035 
130.4324 

134.0639 
137.6980 
141.3343 
144.9729 
148.6135 
152.2560 
155.9003 
159.5^63 



108 
112 
116 
120 
124 
128 
132 
136' 
i4o- 
'144 
148' 
152 
156 
160 
164 
168 
172 
176' 
180 
i84 
188. 
192 
196 

200 
204 



159.5463 

163.1939 
166.8431 

170.4936 

174.1456 

177.7987 
•181.4532 

185.1087 

188.7654 

192.4230 

196.0817 

199.7^13 
203.4018 

207.0631 
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212 

216 

220 

224 

228 

232* 

236 

24o 
244 
248 
252 
256 
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TlovA/ CHAViT FOR J^(x) 







CODING 


SHEETS 




0000.0 


+15-7764.0 
+57.0012.0 


SAX 
TRA 


0020.0 +35.7772.0 
+57.0034.0 


FST 
TRA 


0001 .0 


+00.7774.0 

+22.7772.0 


CLA 
DIV 


0021 .0 -00.0000.0 
-00.0000.0 





0002.0 


+50.0115.0 
+57-0113-0 


TIE 
TRA 


0022.0 -00.0000,0 
-00.0000.0 





0003.0 


-00.0000.0 
-00.0000.0 





0023.0 +00.0000.0 
-00.0000.0 





0004.0 


+07.7772.0 
+57.0172.0 


FMP 
TRA 


0024.0 +35.7766.0 
+57.0040.0' 


FST 
TRA 


0005.0 


+30.0010.0 
+57-0016.0 


FCA 
TRA 


0025.0 +30.7776.0 
+57.01 42.0 


FCA 
TRA 


0006.0 


-00.0000.0 
-00.0000,0 


-«. 


0026.0 +06.7774.0 
+57.0055.0 


FSB 
TRA 


0007.0 


+60.7770.0 
+57.0C25.0 


STO 
TRA 


0027.0 +00.0000.0 
-00.0000.1 


— — — 


0010.0 


+00.0000.0 
-00.0000.0 





0030.0 +30.7776.0 
+57.oo4i.o 


FCA 
TRA 


0011 .0 


+00.0000.0 
-(X). 0000.0 





0031.0 +30.7760.0 
+57.0144.0 


FCA 
TRA 


0012.0 


+01 .0015.0 

+57.0126.0 


ADD 
TRA 


0032.0 +30.7776.0 
+57.0045,0* 


FCA 
TRA 


0013-0 


+00.0000.0 
-00.0000.0 





0033.0 +35.7760.0 

• +57.7772.1 


FST 
TRA 


0014.0 


+05.7776.0 
+57-0103.0 


FDV 
TRA 


0034.0 +30. ([66.0 
+57. 0044.0' 


FCA 
TRA 


0015.0 


+30.7760.0 
+57-0001 .1 


FCA 
TRA 


0035.0 +03.0037.0 
• +57.0146.0 


SUB 
TRA 


0016.0 


+43.0000.0 
+57.0225.0 


XAR 
TRA 


0036.0 '+43.0000.0 
+57.0046.0 


XAR 
TRA 


0017.0 


+35.7762.0 
+57.0031.0 


FST 
TRA 


0037.0 +57.0001 .0 
-^0.7535.0 


..... 
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00^0.0 ■*• 30. 7776.0 rcA 0060.04-43.0000.0 xar 

+57.0050.0 TRA +57.0001.0 TRA 

00^11.0+35.777^.0 FST 0061.0+36.7773.1 DIS 

. +57.0052.0 TRA '+57.0221.0 TRA 

OOte.O +06.777^.0 FSB 0062.0+30.776^1.0 FCA 

+57.CX)71.0 TRA +57.0672.0 TRA 

00^13.0 +35.0010.0 FST 0063.0 +30.7772.0 FCA 

+57.0057.0 TRA +57.0074.0 TRA 

0044.0 +35.7776.0 FST 0064.0 +35.7776.0 FST 

+57.CX)54.0 TRA '+57. 01 74.0* TRA 

0045.0+35.7774.0 FST 0065.0+06.7760.0 FSB 

+57.00^2.0 TRA +57.0153.0 TRA 

0046.0+66.0250.0 CTv 0066.0 +o4. 001 0.0 fad 

+57.0065.0 TRA +57.0143.0 TRA 

0047.0+03.0051.0 SUB 0067.0+35.7764.0 FST 

+57. 01 60.0 TRA +57.0102.0 TRA 



0050.0+35.7774.0 FST 0070.0+00.0004.0 -— 

+57.0063.0 TRA -00.0010.0 -— 

0051 .0 +30.7607. 1 «— 0071 .0 +35.7764.0 FST 

+77.7777.1 -— +57.0105.0 TRA 

0052.0+30.7762.0 FCA 0072.0+35.7776.0 FST 

+57.0064.0 TRA +57.6004.0* TRA 

0053.0+60.7773.0 STO 0073.0+06.0076.0 FSB 

+57.0001.0 TRA +57.0120.0- TRA 

0054.0 +07.7772.0 FMP 0074.0 +06.0076.0 FSB 

+57.0042.0 TRA +57.0020.0- TRA 

0055.0 +35.7760.0 FST 0075.0 +00.0000.0 — 

+57.0066.0 TRA -00.0000.0* • 

CK>56.0 +(X). 7770.0 CLA 0076.0 +00.0000.0 -— 

+57.0166.0 TRA --00.0000.0 — 

CO57.O +30.7772.0 FCA 0077.0 +00.0000.0 — 

+57.0073.0 TRA -00.0000.0 — 



0100.0 
0101 .0 
0102.0 

0103.0 
oioH.o 
0105.0 
0106.0 
0107.0 



+30.7772. 

+57.0110. 

-00.0000. 
-00.0000. 

+30.7766. 

+57.001 il. 

+35.7766. 

+57.011^1. 

+0H.7772. 
+57.0130. 

+0^^.0010. 

+57.00143. 

+i|0.0000. 
-00.0000. 

+00.0000. 
-00.002^. 



FCA 

TRA 

— 

— - 

FCA 

TRA 

rsT 

TRA 

FAD 

TRA 

FAD 

TRA 

- — 

- — 

— 

— 



0120.0 +35.7772.0 FST 

+57.0032.0 TRA 

0121 .0 -00.0000.0 — 

-00.0000.0 - — 

0122.0 +04.0226.0 FAD 

• +57.0171.0 TRA 

0123.0 +35.7762.0 FST 

' +^^J.0^kk.0 TRA 

012^1.0 +41.0001.0 ALS 

+57.0240.0 TRA 

0125.0 +05.7776.0 FDV 

+57.0017.0 TRA 

0126.0 +60.0075.0 STO 

+57.0047.0 TRA 

0127.0 -00.0000.0 — 

-00.0000.0 — 



0110.0 +07.7776.0 FMP 

+57.0176.0 TRA 

0111.0 -00.0000.0 — - 

-00.0000.0 — 

0112.0 +00.7767.0 CLA 

+57.0154.0 TRA 

0113.0 +43.0000.0 XAR 

+57.0131.0 TRA 

0114.0 +30.7762.0 FCA 

+57.0125.0 TRA 

0115.0 +43.0000.0 XAR 

+4i. 0043.0 ALS 

0116.0 +01.7773.0 ADD 

+57.0134.0 TRA 

0117.0 +00.7766.0 CLA 

+57.0135.0 TRA 



0130.0 
0131.0 
0132.0 
0133.0 
0134.0 
0135.0 
0136.0 
0137.0 



+35.7772. 

+57.0140. 

+41.0043. 
+01.7773. 

+4o.ooo4. 
+57.0053. 

+07.0136. 
+57.0122. 

+60.7773. 

+57. 0061. 

+03.7777. 
'+57.0247. 

•+45.5136. 
-65.1770. 

+00.0000. 
-00.0000. 



FST 

TRA 

ALS 

ADD 

ARS 

TRA 

FMP 

TRA 

STO 

TRA 

SUB 

TRA 

— 

— 

— 
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0164.0 +35* 7776.0 fst 

+57.O27I+.O TRA 



01 ^40.0 +07.7776.0 FMP 0160.0 +60.002-^.0 STO 

+57.0026.0 TRA +57.0175.0 TRA 

01^11.0-00.0000.0 —- 0161.0+35.7763.0 FST 

-00.0000.0 -.— -+57.0173.0 TRA 

OIH2.O +35.7774.0 FST 0162.0 +00.7764.0 CLA 

+57.0152.0 TRA ' +57.0145.0 TRA 

OIH3.O +35.0010.0 FST 01 63. 0'+04. 0270.0- FAD 

+57.0056.0 TRA ' '+57.01 75.0 TRA 

0144.0 +05.7776.0 FDV 

+57.0033.0 TRA 

0145.0+35.7760.0 FST 0165.0+30.7764.0 FCA 

+57.0156.0 TRA +57.0200.0 TRA 

0146.0+60.0013.0 STO 0166.0+03.0070.0 SUB 

+57.0162.0 TRA +57.0177.0 TRA 

0147.0+01.7770.0 ADD 0167.0+50.0075.0 TZE 

+57.7752.0 TRA +57-0133.0 TRA 



0150.0+00.0000.0 0170.0+41.0024.0 ALS 

-00.0000.0 — +57.0035.0 TRA 

01^' 1.0 -00.0000.0 — 0171.0+07.7760.0 FMP 

-00.0000.0 — +57.0161.0 TRA 

01^2.0+30.7760.0 FCA 0172.0+06.7774.0 FSB 

+57.0164.0 TRA +57.0^20.0 TRA 

0153.0+35.7776.0 FST 0173.0 +44.7763.0- FSQ 

+57.0165.0 TRA +57.0163.0- TRA 

0154.0+66.0260.0 CTV 0174.0+34.0076.0 FCS 

+57.777' -0 . TRA •• • -+57. 0104.0- TRA 

0155.0+33.0157.0 EXT 0175.0 +04.0106.0 FAD 

+57.0273.0- TRA +57.0272.0 TRA 

0156.0+51.0075.0 TMI 0176.0+06.7774.0- FSB 

+57.0167.0 TRA +57.0024.0- TRA 

0157.0+00.0000.0 — 0177.0+51.0005.0 TMI 

-77.7777.0 — +57.0007.0 TRA 
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0200.0 
0201 .0 
0202.0 
0203.0 
020^1.0 
020^.0 
020^.0 
0207.0 



+05.777^. 
+57.0067. 

+00.0000. 

+6^1.0300. 

+57.0165. 
+55.0300. 

+60.77^6. 
+57.0013. 

+oii.7772. 
+57.0230. 

+60.0150. 
+57.0117. 

+66.0210. 
+57.012I4. 

+00.7766. 
+57.022^. 



FDV 

TRA 

CLA 

CTL 

TRA 

cri 

STO 

TRA 

FAD 

TRA 

STO 

TRA 

CTV 

TRA 

CLA 

TRA 



0220.0 +35.7762.0 FST 

+57.0030.0 TRA 

0221.0 +30.7776.0 FCA 

+57.0241.0 TRA 

0222.0 +64.0230.0 CTL 

+57.0147.0' TRA 

0223.0 -00.0000.0 — 

-00.0000.0 — 

0224.0 +40.0002.0 ARS 

+57.0246.0 TRA 

0225.0 +01.0027.0 ADD 

+57. 0036. 0' TRA 

0226.0 +40.6044.1 -— 

+35.1361.1 — . 

0227.0 +00.0000.0 -— 

-00.0003.1 — 



0210.0 +00.0000.0 

+65.0310.0 

0211.0 +57.0100.0 
+77.0000.0 

0212.0 +00.0000.0 
-00.0005.0 

0213.0 +03.6000.0 
-00.0000.0 

0214.0 +00.0000.0 
-00.0000.0 

0215.0 -00.0000.0 
-00.0000.0 

0216.0 +40.0000.0 
-00.0000.0 

0217.0 +00.0000.0 
-00.0001.0 



CFL 

TRA 
HTR 



0230.0 +35.7772.0 FST 

+57.7770.1 TRA 

0231.0 +35.0076.0 FST 

• +57.0245.0 TRA 

0232.0 +60.7770.0 STO 

+00.7775.0 CLA 

0233.0 +4l .0023.0 ALS 

+60.7760.0- STO 

0234.0 +01.7770.0 ADO 

• *• +60.7770.0 STO 

0235.0 +30.0201 .0 FCA 

+01.7760.0- ADO 

0236.0 +35.0251.0 FST 

+30.7776.0 FCA 

0237.0 +35.0010.0 FST 

+57.0100.0 TRA 



02^40.0 


+60.7775.0 
+57.0256.0 


STO 
TRA. 


0?60.0 +00.0n00.0 
-OC.OOOil.O 


— 


02^11 .0 


+05.7760.0 

+57.0231 .0 


FDV 
TRA 


0261 .0 +40.002i|.0 

+60. [(S^.o- 


ARS 
STO 


02te.O 


+ 35.^/72.0 
+57.0276.0 


FST 

TRA 


0262.0 +03.7766.0 
+51 .0075.0- 


SUB 
TMI 


O2J13.O 


+60.7767.0 
+57.0155.0 


STO 

TRA 


0263.0 +00.7770.0 
+01.7762.0 


CLA 
ADO 


02kk.O 


+00.7766.0 
+57.0257.0 


CLA 
TRA 


0264.0 +60.7767.0 
+03.7766.0 


STO 
SUB 


O2H5.O 


+07.7762.0 
+57.02^2.0 


FMP 
TRA 


0265.0 +51.0207.0 
+00.7766.0 


TMI 
CLA 


O2J16.O 


+ ^5.0000.0 
+57.0123.6 


FNM 
TRA 


0266.0 +01 .7770.0 
+60.7766.0 


ADD 
STO 


O2U7.O 


+U0. 0001.0 

+57.0060.0 


ARS 
TRA 


0267.0 +00.7763.0 

+57.7772.0 


CLA 

TRA 


0250.0 


+00.0003.1 

+77.7770.0 




0270.0 +70.0000.0 

-00.0000.0 


— - 


0251.0 


+00.0000.0 
-00.0000.0 





0271 .0 +00.0000.0 
-00.0001.1 


-— 


0252.0 


+00.0000.0 

-00.0000.0 


... 


0272.0 +if1 .0002.0 

+57.0203.0 


ALS 
TRA 


0253.0 


+00.7771.0 
+03.7770.0 


CLA 
SUB 


0273.0 +60.7762.0 
+57.0112.0 


STO 
TRA 


025^,0 


+51.0023.0 
+42.7772.1 


TMI 
STA 


027^1.0 +3^1.0076.0 

+57. 020^1.0 


PCS 

TRA 


0255.0 


+60.7771.0 

+57.7771.1 


STO 
TRA 


0275.0 +00.0000.0 

-00,0306.0 





0256.0 


+01.0275.0 
+57.0205.0 


ADD 
TRA 


0276.0 +00.7775.0 
+57.0222.0 


CLA 
TRA 


0257.0 


+03.7770.0 

+57.0206.0 


SUB 
TRA 


0277.0 -00.0000.0 
-00.0000.0 


... 



6. 



